This invention relates to a method of operation of a distributed network of computers or a parallel computer, having numerous processors, for determining the evolution ofthe state variables associated with the components of a large physical system, wherein this evolution is describable by a large set of ordinary differential equations, the number of which is equal to the number of state variables, and wherein these equations are subdivided into groups of equations. A vector representing all of the state variables is subdivided into subvectors or blocks, each of which is associated with one of the groups of equations and contains a specified number of vector components. Thereupon, in each group of equations, the components of the associated subvector are treated as unknown functions of time and the remainder of the state variables are lumped into a complementary vector whose components are given estimated values, preparatory to the parallel updating of the state variables with respect to time.
One type of physical system of great interest is an electrical circuit, consisting of many interconnected components such as resistors, capacitors, and inductors. By way of example, Very Large Scale Integrated (VLSI) circuits on a computer chip are representative of such systems. For a judiciously chosen combination of voltages and currents, the changes in these circuit state variables are describable via a system of Ordinary Differential Equations (ODEs). These differential equations express, by way of example, conservation of current in the various circuit nodes, and may include the relationship between currents flowing through inductors and the voltage drops across the inductors in question. In a VLSI circuit, the number of circuit components is very large, giving rise to a similarly large number of differential equations which are to be solved for the unknown voltages and currents associated with the various circuit components. In order to reduce the tremendous amount of computation time, necessitated by the solution of such large numbers of differential equations, it is imperative to solve these equations by parallel computation.
It is noted that the description of the forgoing VLSI circuit can also be used to represent, in terms of mathematically equivalent elements, the components of a mechanical system such as those of a suspension bridge. Indeed, analyses of such bridges are often conducted by way of simulation whereby mechanical structures are represented by electrical components, with currents and voltages representing velocities and forces. In this sense, the set of differential equations may represent the response of either electrical or mechanical systems. It is known from the field of fiuidics that fluidic amplifiers, involving fluid flowing in channels of varying cross sections, can also be described by a system of differential equations of the Forgoing type and that such a fluid system can be simulated by electrical circuit components. In corresponding fashion, as is well known from the simulation of neural networks, a system of biochemical elements is describable by an electrical circuit and the associated set of differential equations. In this sense, differential equations describe the response of fluidic and biochemical systems as well, in a fashion analoguous to that of electrical and mechanical systems.
To explain how the response of a circuit can be described by differential equations, consider the small model of an RLC circuit shown in FIG. 1. The five nodes 1-5 are connected to ground through the five capacitors C.sub.1 -C.sub.5. Three inductors, L.sub.1 -L.sub.3, are in series with two resistors, R.sub.1, R.sub.2, as shown. The directions of the current flows are indicated by arrows. If we number the five nodes by i or j, running from 1 to 5, then the following basic relationships between currents and voltages are valid:
a) A current I, flowing through a resistor of resistance R From node No. i to node No. j, satisfies Ohm's law: EQU I-G(V.sub.i -V.sub.j)=0, (1)
where G=1/R is the admittance of the resistor and V.sub.i, V.sub.j denote the voltages of nodes i and j relative to ground (V=0), respectively.
b) A current I, flowing through a capacitor of capacitance C from node i to node j, satisfies EQU I-Cd/dt(V.sub.i -V.sub.j)=0, (2)
where d/dt denotes the derivative with respect to time.
c) Finally, a current I, flowing across an inductor of inductance L from node i to node j, satisfies EQU V.sub.i -V.sub.j -LdI/dt=0. (3)
The response of the model circuit can be described by the time dependence of the eight state variables V.sub.i, where i runs from 1 to 5, I.sub.L1, I.sub.L2, and I.sub.L3. The eight unknowns satisfy eight differential equations derived by using (1)-(3). The first five of these express the fact that the sum of all currents flowing out of node i minus the sum of all currents flowing into that node, i=1, . . . , 5, respectively, is zero. For example, two currents flow out of node 2, namely the current I.sub.R1 through R.sub.1 into node 3 and the current I.sub.C2 through C.sub.2 to ground, and one current, I.sub.L1, flows from node 1 through L.sub.1 into node 2. Thus, EQU I.sub.R1 +I.sub.C2 -I.sub.L1 =0.
Using (1) and (2), this can be expressed as EQU G.sub.1 (V.sub.2 -V.sub.3)+C.sub.2 V.sub.2 -I.sub.L1 =0,
where V.sub.2 is a simplified notation for the time derivative of V.sub.2 and where G.sub.1 =C.sub.2 =2.
The last three differential equations represent the so-called branch relations (3) associated with the three inductors L.sub.1,L.sub.4, and L.sub.3, respectively. The entire system of differential equations can be written in matrix--vector form: EQU DY-AY=0. (4)
Here, Y is a vector whose components are the eight unknown state variables, D is a diagonal matrix whose diagonal elements are (1,2,3,1,1,-1,-1,-1), and the coefficient matrix ##EQU1## contains the information defining the interaction between the circuit state variables governing the differential equations. The numerical values of the elements of D and A are determined from the admittances, capacitances, and inductances of the circuit of FIG. 1. The example of the current balance at node 2 considered above represents the second equation of the system (4). It should be noted that if in a circuit an inductor is replaced by a resistor, then the numbers of differential equations and unknown state variables diminish by one each.
In state-of-the-art circuit analysis and design, the systems of differential equations to be solved are very large. Furthermore, if the circuit response is to be optimized over choices of parameters (e.g. the values of R, L, and C) those systems must be solved repeatedly, once for each set of parameters.
Very often, the differential equations associated with circuits are stiff, i.e. they have both large and small time constants. Mathematically, for the model problem (4), the time constants are defined as the negative of the reciprocals of the real parts of the eigenvalues of the matrix D.sup.-1 A of equation (4). The model problem becomes stiff if, for example, C.sub.5 =1 is replaced by C.sub.5 =0.01. Efficient numerical integration methods for stiff circuit problems are necessarily implicit, coupling all unknowns in updating the solution from one time level to the next. The updating thus requires, at each time point, the solution of very large systems of linear or nonlinear algebraic equations.
The computational task outlined above cannot be accomplished fast enough on a conventional serial computer. It is therefore imperative to make use of parallel processing on a distributed network of computers. However, designing schemes capable of exploiting parallelism at the level of the numerical integration method is difficult because the problems to be solved are intrinsically sequential in time. As a consequence, only limited progress has been made in this endeavor thus far. This is apparent in studying the Prior Art discussed hereafter.